#include "ErrorEstimator.h"
#include <string>
#include <cstring>
#include <sstream>
#include <iostream>
#include <math.h>
#include <cmath>


ErrorEstimator::ErrorEstimator (pMesh theMesh){

}


void ErrorEstimator::LeitorDeLista (ErrorAnalysis *pErrorAnalysis, pMesh theMesh){

	std::list<pEntity> elementList;
	pErrorAnalysis->getRefUnrefElementsList(theMesh,elementList);
	
	set<int> ListofElements;

  
  	int i = 0;
	for (std::list<pEntity>::iterator iter = elementList.begin(); iter!=elementList.end(); iter++){
	//cout << i++ << "\t" << EN_id ( (*iter)->get(0,0) ) << " " << EN_id ( (*iter)->get(0,1) ) << " " << EN_id ( (*iter)->get(0,2) ) << "\t" << std::setw(5) << pErrorAnalysis->getLevelOfRefinement(*iter) << endl;
	ListedFaces.insert(*iter);
	
	EN_attachDataDbl (*iter, MD_lookupMeshDataId( "Ref_factor" ), pErrorAnalysis->getLevelOfRefinement(*iter));

	}
}







void ErrorEstimator::VertexList (pMesh theMesh){ // Pega os elementos do set de faces a remover e levanta todos os vertex

	set<pFace>::iterator itfaces;
	for (itfaces=ListedFaces.begin() ; itfaces != ListedFaces.end(); itfaces++ ){
		pVertex vt;
		for(int y=0;y<=2;y++){
			vt=F_vertex(*itfaces,y);
			AllNodes.insert(vt);
		}
	}

	cout << "A malha possui " << M_numVertices(theMesh) << " nodes ao total" << endl;
	cout << "AllNodes.size(): " << AllNodes.size() << endl;
}









void ErrorEstimator::Calc_Factor (){ 

	set<pFace>::iterator itfacs;
	double NHD = 0;
	double F = 0;
	double D=1;
	for(itfacs=ListedFaces.begin(); itfacs!=ListedFaces.end(); itfacs++){
		EN_getDataDbl (*itfacs, MD_lookupMeshDataId( "Ref_factor" ), &NHD);

		
		D = pow(2,NHD);
		//cout << endl << "O denominador é: " << D << endl << endl;
		
		F = 1/D;
		//cout << endl << "O fator calculado é: " << F << endl << endl;
		
		EN_attachDataDbl (*itfacs, MD_lookupMeshDataId( "Size_factor" ), F);
	}

}


void ErrorEstimator::Calc_Cl (pMesh theMesh, double Factor){  

	float h;
	float cl;
	int contadordeCLS=0;

	set<pFace> Vfaces;
	set<pEdge> Vedges;
	set<pFace>::iterator itfacs;
	set<pVertex>::iterator itVertx;
	set<pEdge>::iterator iteds;

	
	
	
	for(itVertx=AllNodes.begin(); itVertx!=AllNodes.end(); itVertx++){ // Percorrendo o set AllNodes e levantando todas as faces vizinhas	

	

		int Numfaces = V_numFaces(*itVertx);
		int Numedges = V_numEdges(*itVertx);


		for(int y=0;y<Numedges;y++){
			pEdge ted = V_edge(*itVertx,y);
			Vedges.insert(ted);
		}


		for(iteds=Vedges.begin();iteds!=Vedges.end(); iteds++){
			int Efaces = E_numFaces(*iteds);
			for(int e=0; e<Efaces; e++){
				pFace tfc = E_face(*iteds, e);
				Vfaces.insert(tfc);
			}

		}


		double h=0;

		for(itfacs=Vfaces.begin();itfacs!=Vfaces.end(); itfacs++){

			double SF = 0;
			EN_getDataDbl (*itfacs, MD_lookupMeshDataId( "Size_factor" ), &SF); // Pegando o Fator de tamanho do triângulo

			// Medindo o tamanho do H do ponto

			pVertex pts[2];
			pVertex vertx;
			int V = 0;
			set<int> pontos;

			for (int i=0; i<=2; i++){ // levantando todos os pontos...
				vertx = F_vertex(*itfacs, i);
				if (EN_id(vertx)!=EN_id(*itVertx)){
					pts[V]=vertx;
					V++;
				}
			}

			double xyz1[3];
			double xyz2[3];
			double xyz3[3];

			double x0 = xyz3[0];
			double y0 = xyz3[1];

			V_coord(pts[0], xyz1);
			V_coord(pts[1], xyz2);
			V_coord(*itVertx, xyz3);

			double m = (xyz2[1]-xyz1[1])/(xyz2[0]-xyz1[0]);
			double c = -m*xyz2[0]+xyz2[1];
			double pc = m*x0-y0+c;
			pc = sqrt(pc*pc);
			double k = pc/sqrt((m*m)+1);
			k = sqrt(k*k);
			h=h+(k*SF); // Somando os h's multiplicados pelos fatores de tamanho


		}

		Vedges.clear();
		Vfaces.clear();
		cl = h/Numfaces;
		VCL.insert(std::pair<int,float>(EN_id(*itVertx),cl)); // Inserindo o CL no map VCL
		//cout << endl << "INSERIDO CL NO MAP --- " << cl << endl; 
		contadordeCLS++;
	}
	cout << endl << "O map tem " << contadordeCLS << " CLS" << endl << endl;
}



void ErrorEstimator::SaveFileListCls (){


	system("rm ListaCls");
	ofstream Myfile("./ListaCls", ios::app);
	
	if (!Myfile) {
		Myfile.open("./ListaCls");
	}

	// ROTINA QUE CRIA A LISTA DE PONTOS INICIAL DO ARQUIVO - IDS E COORDENADAS DE TODOS OS PONTOS

	Myfile << "$CLS"<< endl;

	map<int,float>::iterator itVCls;

	for(itVCls=VCL.begin();itVCls!=VCL.end();itVCls++){ // lendo o map e gravando no arquivo os dados.

		Myfile<<(*itVCls).second<<endl;

	}

	Myfile << "$NODES"<< endl;

	for(itVCls=VCL.begin();itVCls!=VCL.end();itVCls++){ // lendo o map e gravando no arquivo os dados.

		Myfile<<(*itVCls).first<<endl;

	}

	Myfile.close();

}

void ErrorEstimator::SaveFileListFaces (){

	system("rm ListadeElementos");
	ofstream Myfile("./ListadeElementos", ios::app);
	
	if (!Myfile) {
		Myfile.open("./ListadeElementos");
	}



	set<pFace>::iterator itfacs;

	for (itfacs=ListedFaces.begin() ; itfacs != ListedFaces.end(); itfacs++ ){ // Esse for aqui percorre o set de faces

		Myfile << EN_id(*itfacs) << endl;

	}

	Myfile.close();

}


















